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Abstract 

Quantum molecular dynamics is applied to study the ground state properties 
of nuclear matter at subsaturation densities. Clustering effects are observed 
as to soften the equation of state at these densities. The structure of nuclear 
matter at subsaturation density shows some exotic shapes with variation of 
the density. 
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I. INTRODUCTION 



One of the main interests of heavy-ion physics and astrophysics is the property of nuclear 
matter in extreme conditions. Its high-density behavior is important for the scenario of 
supernova explosions, the evolution of neutron stars, the reaction process of high-energy 
heavy-ion collisions, quark-gluon plasma and so on. Properties of the nuclear matter below 
saturation density, on the other hand, are essential in describing the multi-fragmentation in 
the heavy-ion collisions, the collapsing stages in supernova explosions and the structure of 
neutron star crusts. Here the saturation density is the density of the energy-minimum state 
of the nuclear matter at a fixed proton ratio. For symmetric nuclear matter, the saturation 
density is the normal nuclear density po = 0.165 fm^^. Since the matter is unstable below 
the saturation density, an inhomogeneous state is expected to appear below the saturation 
density. 

Besides, supernova matter (SNM) and neutron star matter (NSM) are also interesting 
from the viewpoint of the nuclear shape. At low densities, nuclei in these matters are 
expected to be crystalized so as to minimize the long range Coulomb energies. They melt 
into uniform matter at a certain density close to the saturation density. Then what happens 
in between? About a decade ago, three groups ||T]-^ suggested that nuclei have exotic 
structures in the SNM and/or NSM. They showed that the stable nuclear shape changes from 
sphere to cylinder, slab, cylindrical hole, and to spherical hole with increase of the matter 
density. The favorable nuclear shape is determined by a balance between the surface and 
Coulomb energies. In the liquid drop model, a simple geometrical argument demonstrates 
that the favorable shape changes as a function of the volume fraction of the nucleus in the 
cell, independently of specific nuclear interactions [HQ. /,From recent studies in the liquid 
drop model and the Thomas Fermi calculations [^0, the non-spherical shapes of nuclei 
are expected in a density range at about half the saturation density although the range 
depends on the choice of nuclear interaction. Furthermore these shapes are expected to 
survive even if the shell effects are taken into account 
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These exotic nuclear shapes may cause substantial astrophysical consequences. In SNM, 
neutrino absorption by nuclei may modify leptonic energy of the matter and lead to a 
significant change in the strength of supernova explosions ||^. In neutron stars, the exotic 
nuclei may also affect the pinning of superfluid neutron vortices to nuclei in inner crusts of 
neutron stars, which are considered to be the initial step of neutron star glitches. Recently, 
Mochizuki and Izuyama demonstrated that the exotic nuclear shapes actually play an 
essential role in trapping vortices in a microscopic mechanism of the glitches. 

The spatial fluctuation is important in describing properties of SNM and NSM because 
geometrical distributions of nucleons may affect neutrino reaction rate in SNM and inter- 
actions with vortices in NSM. However, there have been only two works P,p!0| which allow 



arbitrary nuclear shapes in their models. These two works 0,|TO| use the Thomas Fermi 



approximation and treat SNM while none has been done for NSM. Williams and Koonin 
have investigated the structure of neutron-proton symmetric nuclear matter with the proton 



ratio Zj A = 0.5, while the authors of Ref. |T0[ studied the matter with Z/A = 0.285 at the 



total entropy per baryon s = 1. Unfortunately, some spatial fluctuation due to a possible 
cluster correlation or non uniform distribution is neglected since their calculations with the 
Thomas-Fermi approximation are based on a one-body treatment of the matter. 

In the field of nuclear reaction study, molecular dynamics has become one of the most 
powerful approach to simulate the fragmentation. Advantages of molecular dynamics for the 
investigation of heavy-ion reactions are that no reaction mechanism is required to assume 
and that the fluctuation of the system is automatically included. 

Within various kinds of molecular dynamics, the quantum molecular dynamics (QMD) 
TTHl6[ approach has been proposed to study high energy heavy-ion collisions. QMD has 



been also used for the analysis of fusion reaction, nucleon-induced reaction, fragmentation 
in collisions between heavy systems, and so on. 

In this paper, we apply the QMD method to the investigation of the equation of state 
(EOS) and the structure of nuclear matter at subsaturation densities. A similar calculation 



using QMD is done in Ref. aiming to see the EOS of the nuclear matter. It was 



reported that the clustering in the matter significantly softens the EOS below saturation 
density. They have used, however, only 254 nucleons in a cell to simulate finite-temperature 
nuclear matter. This number seems insufficient to investigate the nuclear matter properties. 
We then perform the QMD calculation by using much more particles in a cell, and see the 
EOS of nuclear matter at zero-temperature. 

In the present QMD, we introduce the Pauli potential and the momentum-dependent 
interaction in order to simulate the Fermionic feature in a phenomenological way and the 
energy-dependence of the optical potential. The binding energies of finite nuclei and the 
saturation properties of the nuclear matter are well adjusted to the empirical values. As 
a result of this new version of the QMD, we can approximately describe the properties of 
neutron-rich isotopes such as the density profile of ^^Li and the effects of surface neutrons in 



the low-energy collisions . With these ingredients of the QMD and the periodic boundary 
condition, we investigate the ground state properties of the nuclear matter at subsaturation 
densities. 

In Sect. we describe our model based on the QMD. In Sect. |ITT|, the EOS and the 
structure of nuclear matter is discussed. Finally, summary is given in Sect. |V[ 

II. QMD SIMULATION FOR NUCLEAR MATTER 



The present QMD code is an improved version of our previous one |]I6[. In the previ- 
ous works |jl6|,|l8|, where we have employed very simple effective interactions in the QMD, 
we have nicely reproduced several observables in the nucleon-induced reactions at energy 
region of about 100 MeV~5 GeV. Several improvements of this QMD are done to enable 
the calculations for the ground state of the nuclear matter for wide range of density and 
proton ratio. First we have introduced a phenomenological potential, namely Pauli poten- 
tial to describe the Fermionic feature of nucleons. By this Pauli potential, we can uniquely 
determine the ground state of the finite nuclei and the nuclear matter by seeking the energy 
minimum state. Second we have introduced momentum-dependent interaction, which is also 
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an important feature of the Fermion system with finite range interaction. In the following, 
we explain the detail of the present QMD and how it describes the nuclear matter. 



A. The total wave-function and the equation of motion 



In QMD, each nucleon state is represented by a Gaussian wave-function of width L, 



(27rL)3/4 



exp 



(r - R.)^ 
4L 



n 



(1) 



where and Pj are the centers of position and momentum of i-th nucleon, respectively. 
The total wave-function is assumed to be a direct product of these wave-functions. Thus the 
one-body distribution function is obtained by the Wigner transform of the wave-function. 



/(r,p) = I3/i(r,p), 



(2) 



/i(r,p) = 8-exp 



(r-R,)2 2L(p-F,f 



2L 

The equation of motion of Rj and Pj is given by the Newtonian equation 

dH . dH 



(3) 



R, 



(4) 



and the stochastic N-N collision term. Hamiltonian H consists of the kinetic energy and the 
energy of the two-body effective interactions. 



B. Effective interactions 

The Hamiltonian is separated into several parts as follows, 

H^T + Vpauli + Viocal + VmD , (5) 

where T, Vpauii, Mocai and Vmd are the kinetic energy, the Pauli potential, the local 
(momentum-independent) potential and the momentum-dependent potential parts, respec- 
tively. 

5 



The Pauli potential p^ . |T9| -pT| is introduced for the sake of simulating Fermionic proper- 
ties in a semiclassical way. This phenomenological potential prohibits nucleons of the same 
spin a and isospin r from coming close to each other in the phase space. Here we employ 
the Gaussian form of the Pauli potential as 



V -'c ( 

2 \qoPo, 



(6) 



2g2 2pl 

In the local potential part we adopt the Skyrme type with the Coulomb and the symmetry 
terms as explained in Eq. (5) of Ref. |T 



local 



'oY (l + r)p5^ 



+ ~^Y1 CiCj ff d^Tid^rj- — ^ — -pi(ri)pj(rj) 
C 

+ ^ E (l-2|ci-Cj|) pij. 
2P0 ij(^i) 



(7) 



In the above equation, q is 1 for protons and for neutrons, while < pi > is an overlap of 
density with other nucleons defined as 



< Pi> = J2 Pij = J2 Pi(^) Pj(^) 



(8) 



The momentum-dependent term is introduced as a Fock term of the Yukawa-type inter- 
action. We divide this interaction into two ranges so as to fit the effective mass and the 
energy dependence of the real part of the optical potential, as 



MD 



T^(l) , ta(2) 
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2Po i 1 + 

The parametrization of the constants included in the above effective interactions will be 
discussed in Sec. |11K| . 
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C. Energy minimum state 



With inclusion of the Pauh potential, we can define the ground state as an energy- 
minimum state of the system. To get the energy-minimum configuration, we use the following 
damping equation of motion, 

■ dH OH ■ OH dH . 

where /xr, and /xp are the damping coefficients with negative values when we need to cool 
the system. 

We first distribute the particles randomly in the phase space and cool down the system 
according to the damping equation of motion until the energy reaches the minimum value. 
Sometimes the system stops at the local minimum during the cooling. We thus try again and 
again this cooling procedure with a different initial state and seek the real energy minimum 
state. 

For finite nucleus and infinite system above saturation density, this procedure works 
rather well. For infinite system at subsaturation densities, however, there are many local 
minimum states around the real ground state, which differ from the ground state in the 
detail of the surface configuration of the clusters. Since the energy difference from the 
ground state is the order of 10 keV/nucleon in this case, we accept these states as ground 
states and neglect the small differences of the configuration. 

D. Periodic boundary condition 

In order to simulate infinite nuclear matter with finite number of particles, we use a cubic 
cell with periodic boundary condition. The size of the cell is determined from the average 
density and the particle number. The periodic boundary condition is introduced as follows: 
As shown in Fig. |1|, we prepare 26 (= 3'^ — 1) surrounding cells, where the particle distribution 
refiects the distribution of the central cell exactly. The particles in the central cell move 
according to the interaction with all particles in the same cell and in the surrounding cells 
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as well. The particles in the surrounding cells obey exactly the same motions as those in 
the central cell. Thus the Hamiltonian per cell is written as 



H 



T,+ H'lf (R, - + Dceii, P^, P,) + ■ ■ 



coll=0,---,26 
j = l,---,iV 



(11) 



('21 

where is one-body part (kinetic energy), Hlj is the two body part of the Hamiltonian 
and Dceii are the relative position of surrounding cells from the center. Note that the indices 
"cell" runs from (the center cell) to 26 (surrounding cells) and Dq = 0. 



E. Parametrization of the constants 

We have twelve parameters in the effective interactions of the Hamiltonian Eq. (|^), i.e., 
Cp,go5Po5 T, Cs, C^^ ,C^^ , f^i, and the Gaussian width L. We should parametrize 
these constants to reproduce properties of the ground states of the finite nuclei and saturation 
properties of the nuclear matter. 

We first determine the parameters of Pauli potential qQ,Po and Cp, apart from the other 
effective interactions, by fitting the kinetic energy of the exact Fermi gas at zero temperature 
and at various densities. For this, we define the free Fermi gas system as a ground state for 
the Hamiltonian including only the kinetic energy and the Pauli potential by making use of 
the damping equation of motion Eq. ([T0|) and the periodic boundary condition with 1024 
particles in a cell. In Fig. |^, we show the kinetic energies (the solid squares) and the total 
energies (the open squares) obtained by using a parameter set for the Pauli potential as 

Cp = 207 MeV, po = 120 MeV/c, go = 1-644 fm. (12) 

In the same figure we draw the exact energy of the Fermi gas by the solid line. Although 
there are some other parameter sets which can reproduce the exact energies of the Fermi 
gas in the same form, i.g. that used in Ref. |T^, we choose the above parameter set to 
get good properties of the ground states of the finite nuclei with other effective interactions 
particularly in combination with the momentum-dependent interaction. 
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Among remaining nine conditions, four are attributed to the momentum-dependent in- 
teraction as follows. We calculate the single particle potential of momentum p in ideal 
nuclear matter at the normal nuclear density, which leads to 



U{p,Po) = f/iocai + Umd{p) 

^4 x-i 



1 + 



p p 



+ 



1 + 



p-p 
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X 



X 



■ (14) 



We fit the energy-dependence of this potential to the experimental data. In Fig. ^, we plot 
the energy dependence of the real part of the optical potential (the open circles and squares) 
obtained from the experimental data of Hama et al. p2| for p-nucleus elastic scattering. 
iFrom this figure, we pick up three conditions, i.e. U (0) = —80 MeV, U{p) = at i^iab = 200 
MeV, and U{p ^ oo) = a + (3 = 77 MeV. For another condition, we use the value of effective 
mass defined by 



1 _ 1 ^ IdUuD 
m* m \p dp 



(15) 



P=PF 



We take the value of m* = 0.8 m at p = po- 

Other three conditions are coming from the saturation condition, i.e. the energy per 
nucleon E/A = —16 MeV at p = po (0.165 fm~^) and the value of incompressibility K. 

The last two parameters are given by hand. One is the value of the coefficient of the 
symmetry term Cg. We take 25 MeV for Cs- This value leads to the symmetry energy 34.6 
MeV for the nuclear matter at the saturation density (see Sec. [Ill Bp . The other is the width 
of the Gaussian wave packet L, which is a free parameter in QMD model. This value affects 
ground state properties of finite nuclei and infinite nuclear matter below saturation densities, 
while it does not change those of infinite nuclear matter above saturation densities. We then 
choose this value to give nice fitting the binding energies of finite nuclei. 
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It should be noted here that we cannot determine these parameters from the above 
conditions in an analytical way, since the Fermi distribution is not exactly achieved by 
the Pauli potential and the additional potential energy included in the Pauli potential. In 
addition, the saturation properties of the nuclear matter should be realized in the simulated 
matter for the main purpose of this paper. We then simulate the nuclear matter by the 
QMD with the periodic boundary condition using 1024 particles in a cell. We search the 
energy minimum state by the damping equation of motion as discussed above and adjust the 
parameters. By this method, we have fixed three parameter sets according to the value of 
incompressibility which are shown in Table 1. We have prepared three kinds of equation 
of state, namely Soft (ir=210 MeV), Medium {K=2m MeV) and Hard (ir=380 MeV) EOS. 
These values of incompressibility K are subtracted from the results of EOS (shown in Fig. ^ 
by fitting its curvature at the saturation point to the following parabolic form, 

E/A = (P - Po)' - 16 . (16) 

This parabola is also shown in Fig. |^. 

Here the single particle potential shown in Fig. |^ are also calculated by the simulated 
nuclear matter with Pauli potential and other effective interactions instead of ideal nuclear 
matter. The results are denoted by the crosses in Fig. ^ and well coincident with the results 
of ideal nuclear matter except for the low energy part, where the Pauli potential is effective. 
Though this result in Fig. |^ is obtained with a parameter set of Medium EOS, results with 
Soft and Hard EOS are same as Medium EOS within 2 MeV for all energy region. 

In Fig. ^, we plot the binding energies of the ground state of finite nuclei obtained by the 
damping equation of motion Eq. (|T0|) with three parameter sets, i.e.. Soft (the long dashed 
line). Medium (the dashed line) and Hard (the solid line) EOS. All of them reproduce well 
the global trend of the binding energies of various nuclei except for light nuclei from ^^C 
to ^°Ne. It might be due to the specific structures of these light nuclei, which are not well 
described by the present QMD. 
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F. "Screened" Coulomb potential 



For the neutron star, the same number of electrons exist as proton, since the nuclear 
matter in the neutron star should be charge-neutral. Hence the energy of the system remains 
finite even if we calculate the Coulomb interaction of protons and electrons. However, 
Coulomb interaction has so long range that the Coulomb energy depends on the cell-size in 
our treatment of the infinite system by the periodic boundary condition with the surrounding 
neighbor cells. To avoid this cell-size dependence, we introduce a cutoff of the Coulomb 
interaction in the way of "screened" Coulomb potential. We use the following "screened" 
Coulomb interaction instead of the second term of Eq. (J^) for the nuclear matter calculations, 

Vc' = y E jj ciV.ciV /"P[-'\7f/"-] p.(r.)p,(r,) . (17) 

In this equation, rscr is the "screening" length, which we use 10 fm in the present study. 
The physical screening of the Coulomb potential by the electron localization is, however, 
estimated to be much larger in the case of normal nuclear density . Thus our "screening" 
should be considered as a technical approximation to avoid this cell-size dependence and to 
make the numerical calculation feasible. For this purpose, Tscv should be smaller than the 
cell size. On the other hand, to keep the proper description of finite nuclei, it should be lager 
than the size of nuclei. By the "screened" Coulomb interaction with rscr = 10 fm, however, 
the binding energies are slightly modified particularly in heavy nuclei. We compare the 
binding energies obtained by the "screened" Coulomb interaction with that of the normal 
one in Fig. ^ The dotted line is the result of the "screened" Coulomb interaction with 
Medium EOS. Though the binding energies of heavy nuclei increase, the binding energies 
of finite nuclei still have the maximum around at A ^ 100. This feature is important to 
describe clustering of the matter at low densities. 
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III. EQUATION OF STATE AND THE STRUCTURE OF NUCLEAR MATTER 

In this section we study properties of nuclear matter at several conditions. It is desirable 
to use a cell large enough to include several periods of structure and to avoid the spurious 
effects of boundary condition on the structure of the matter. Though our calculation with 
typically 1024 particles in a cell is not fully satisfactory in this respect, we consider it is 
enough for semiqualitative discussions at the beginning of this study. Actually, the global 
quantities, i.g. the ground state energy of the system, is well saturated at this number 
of particles in a cell. In Fig. |^, we show the energy per nucleon of the infinite system as a 
function of the particle number in a cell for four average densities from p = 0.4po to 2.2po- For 
all densities, the energy of the system have already approached the asymptotic value above 
256 particles within 100 keV/nucleon. In this study we then simulate the infinite system by 
the periodic boundary condition with 1024 or 2048 particles in a cell, and investigate ground 
state properties of the nuclear matter. 

A. Symmetric nuclear matter 

We first perform the calculations for symmetric {Z/A = 0.5) matter at zero temperature 
to simulate supernova matter (SNM) in the collapsing stage. Figure ^ shows the energy per 
nucleon as a function of the average density. The solid squares indicate the energy of "uni- 
form" nuclear matter while open squares are the results of energy-minimum configurations. 

The "uniform" matter energy is calculated as follows: First we distribute nucleons ran- 
domly and cool the system only with the Pauli potential. Pauli potential is repulsive and 
does not spoil the uniformity of the system. Then we impose the other effective interactions 
and cool only in the momentum space fixing the positions of particles. The system turns 
out to be approximately uniform with this procedure. Note that the simulated "uniform" 
matter is not exactly the same as ideal nuclear matter since the latter is continuous and 
completely uniform. 
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Both cases of uniform and energy-minimum configurations have almost the same energy 
per nucleon for the higher densities as is seen in this figure. Below saturation density po, the 
energy per nucleon of the energy-minimum configuration is lower than the uniform case. The 
deviation amounts to about 5 MeV. As we see in the following, this is due to the structure 
change of the matter from uniform to non uniform structure such as clusterized one. 

This change of the structure is displayed in Fig. ^ Above 0.8 po, the system is almost 
uniform and no specific structure is seen. Below 0.8 po, however, there appear some voids 
between the matter. As the density gets lower, the voids develop and nuclei are surrounded 
by the voids. Below 0.2 po each nucleus is separate, while above that density nuclei are 
connected to form some transient structures. This change of the structure is basically the 
same as reported in the previous works Furthermore, the transient shape of the 

matter like hole, slab and cylinder are partly seen in our calculation although the nuclear 
surface shape in QMD is somewhat complicated. It should be noted that the nuclear shapes 
do not show exact symmetry properties (sphere, cylinder or slab) assumed in the previous 
liquid-drop and Thomas- Fermi studies where the assumption of the symmetries 

leads to clear changes of the nuclear shape. Our result suggests that during the transition 
from homogeneous to inhomogeneous matter, the nuclear shape may not have these simple 
symmetry properties. It is possible that this is due to the incomplete minimization of the 
energy because there are several local minimums around the real ground state and the energy 
difference due to the nuclear shape is extremely small. However, these shapes are not so 
strange in the transient SNM at the collapsing stage because the matter is not in perfect 
equilibrium. 

The change of the shape obtained in the QMD calculation may affect the neutrino reac- 
tion rate in SNM. It has been pointed out that the nuclear distributions are essentially 
determined by a delicate balance (of the order of 1 keV/nucleon) between the surface and 
Coulomb energies. In the present treatment, however, we neglect such a tiny energy differ- 
ence between the nuclear shapes around the ground state. We then need further improvement 
of our treatment to investigate the affect of the shape on the neutrino reaction rate in SNM. 
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Nevertheless, our results provide the EOS of the symmetric matter with sufficient accu- 
racy as well as global nuclear structure in the matter. 

B. Asymmetric nuclear matter 

In the neutron star matter (NSM), the beta-equilibrium is achieved and the proton ratio 
is given by the energy-minimum condition. The left panel of Fig. || is the energy per nucleon 
of nuclear matter with several proton ratios. The electron kinetic energy is not included in 
this figure, though the Coulomb interaction of protons and uniform electron background is 
included. The energy per nucleon at p = po is fitted by 

E/A ^ -16.2 + 34.6^ — —- [MeV]. (18) 

In other words, the symmetry energy at normal density is 34.6 MeV in our calculation. The 
symmetry energies at p = 0.6po and 0.2po obtained in the same way are 23.0 and 18.9 MeV, 
respectively. 

Including the electron kinetic energy at zero temperature, we get the total energy of 
NSM as a function of the proton ratio shown in the right panel of Fig. It can be seen 
from this figure that the proton ratio that gives the energy-minimum of the system is 0.032 
(p = 0.2po), 0.043 (p = 0.6po) and 0.072 (p = l.Opo). 

The structure of asymmetric nuclear matter at low density (O.lpo) is shown in Fig. |^. 
Even though the proton ratio is small, nucleons form the cluster structure at low density. 
If the proton ratio is very small, some neutrons can not stay inside the cluster but overflow 
into the space; clusters are floating in the neutron sea. When the proton ratio increases, 
free neutrons are absorbed into the clusters. 

The departure from the spherical symmetry of the nuclear (cluster) shape is also seen 
as in the symmetric matter. This may cause some consequences in the standard scenario of 
pulsar glitches. In the scenario, vortices in the superfluid neutron sea are supposed to be 
pinned to clusters (nuclei) and accumulated in the inner crusts because neutrons are normal 
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in nuclei. The strength of the pinning, which is important in the scenario, depends on the 
geometry of nuclei as well as on the superfluid energy gap. 

Finally, we have to recognize that the simulation of infinite system still needs much lager 
number of particles. At least a cell must include several periods of structure for distinct 
conclusions. If only one or two units of the structure is included in a cell, the size of the 
unit structure is the same or half size of the given cell size. The lattice is also limited to be 
cubic. In this case, the results are dependent on the boundary condition which is artificially 
imposed. However, we emphasize that the QMD framework can be used also for the NSM 
and the symmetric matter as has been used for the reaction studies. Together with some 
refinements of the surface energy, the present method will be able to describe from stable 
and unstable nuclei to homogeneous and inhomogeneous nuclear matter. The problem of 
the computational time in this study is expected to be solved soon. 

IV. SUMMARY 

We have proposed the QMD approach for the description of nuclear matter in wide range 
of density and proton ratio. We can well reproduce the finite nuclear properties for various 
mass range by inclusion of the Pauli potential and the momentum-dependent interaction. 
We have investigated the EOS of nuclear matter by simulating the infinite system with our 
QMD. Below the saturation density, clustering of the system was observed, which softens 
the EOS by lowering the energy per nucleon up to about 5 MeV. 

We have shown the structure of the nuclear matter at subsaturation density. The tran- 
sient shape of the symmetric nuclear matter, such as hole, slab, cylinder and sphere, pre- 
dicted in previous works with analytic model and Thomas- Fermi calculations are 
partially seen in our calculation. However the structure of nuclear matter at subsaturation 
density appears rather vaguely in our case. This result suggests that during the transition 
from homogeneous to inhomogeneous matter, the nuclear shape may not have these sim- 
ple symmetry properties. We need, however, further investigation increasing the particle 
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number to get quantitative conclusion. 

For asymmetric nuclear matter we have obtained the proton ratios ZjA — 0.032 (p = 
0.2/9o), 0.043 (p = 0.6po) and 0.072 (p = l.Opo), which give the energy-minimum of the 
system for the fixed average densities. At considerably low proton ratio, we have observed 
a neutron sea in which the normal nuclei are floating around. 

In this paper we have presented our first results on the infinite nuclear matter by the use 
of molecular dynamics method. Though it is still necessary to enlarge the particle number, 
our results agree quantitatively to previous studies which include much more assumptions 
and restrictions in the models. 

Our model contains further possibility for the simulation of dynamical evolution of infinite 
nuclear matter such as supernova explosion, the glitch of the neutron star and the initial 
stage of the universe. Intensive and systematic study of the nuclear matter with present 
model will be important since it contains less assumption than the foregoing models as to 
the structure of the matter. 
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TABLES 



TABLE I. Effective interaction parameter set 





Soft (a =210 MeV) 


TV fl" J * / oor\ '\ X '\T\ 

Medium (A =280 MeVj 


TT J / oor\ li X '\T\ 

Hard (iv =380 MeV) 


a. (MeV) 


—223.56 


—92.86 


—21.21 


p (MeV) 


298.78 


169.28 


97.93 


r 


I.IdddT 


1.33333 


1.6d6d7 


C<= fMeV) 


25.0 


25.0 


25.0 


dx^ (MeV) 


-258.54 


-258.54 


-258.54 


dx^ (MeV) 


375.6 


375.6 


375.6 


/ii (MeV) 


2.35 


2.35 


2.35 


/X2 (MeV) 


0.4 


0.4 


0.4 


L (fm2) 


2.1 


2.1 


2.05 
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FIGURES 

FIG. 1. An illustrative explanation of cell configuration. 26 surrounding cells (only 8 cells 
are displayed in this figure) have exactly the same distribution of particles as the central cell. The 
relative position vector of each surrounding cell from the central cell is Dceii and Do = 0. 

FIG. 2. Energy per particle of free Fermi gas. The solid line shows the exact value. The cases 
of the molecular dynamics calculation with only Pauli potential are shown by the solid squares 
(kinetic energy) and the open squares (total energy). 

FIG. 3. The energy dependence of the real part of the optical potential. The open circles and 



squares indicate the results obtained from the experimental data of Hama et al. |22] for p- nucleus 
elastic scattering. The solid line denotes the single particle potential calculated by Eq. (|T3|) in 
ideal nuclear matter with the parameter set of Medium EOS. The short and long dashed lines 



show a + f3 + Cexg{x,y) and C^g{x,y) of Eq. (|13D, respectively. The single particle potential 



calculated by the simulated nuclear matter with Pauli potential is shown by the crosses. 

FIG. 4. Binding energies of finite system obtained by the damping equation of motion with 
three parameter sets, i.e., Soft (the long dashed line), Medium (the dashed line) and Hard (the 
solid line) EOS. The solid squares denote experimental data. The dotted line indicates the binding 
energy per nucleon obtained by using the "screened" Coulomb interaction in the case of Medium 
EOS. 

FIG. 5. Particle number dependence of the energy per nucleon of the infinite system for four 
average densities from p = 0.4po to 2.2po- 
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FIG. 6. The energy per nucleon of symmetric nuclear matter (ZjA = 0.5) at zero temperature 
as a function of the average densities. ^Prom the left, the open squares are the results with Soft 
{K = 210 MeV), Medium {K = 280 MeV) and Hard EOS (K = 380 MeV) obtained by the 
damping equation of motion searching the energy-minimum configuration in the full phase-space. 
The solid squares indicate results obtained from the spatially uniform distribution. The kinetic 
energy of the electron is not included. We use 1024 particles in a cell for all cases. 

FIG. 7. The structure of symmetric nuclear matter. ^From the upper left, the average density 
is p = l.Opoj O.SpO) 0.6/90) 0.4/90, 0.2/9o and O.lpo. The white circles denote neutrons and red circles 
are protons. Nuclear potential of Medium EOS is used. We use 2048 particles in a cell for these 
cases and the size of a cell is indicated in the figure. 

FIG. 8. The left panel: the energy per nucleon of nuclear matter as a function of the proton 
ratio Z/A for the three fixed average densities p = 0.2po (the solid triangles), 0.6po (the solid 
squares), and l.Opo (the solid circles). The right panel: same as the left panel but with the kinetic 
energy of the electrons. 

FIG. 9. The structure of asymmetric nuclear matter. ^Prom the upper left, the proton ratio 
Z/A is 0.098, 0.195, 0.293 and 0.391, while the average density is O.lpo for all panels. The white 
circles denote neutrons and red circles are protons. Nuclear potential of Medium EOS is used. We 
use 1024 particles in a cell for these cases and the size of a cell is indicated in the figure. 
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